import Jama.Matrix;
import java.util.ArrayList;
import java.util.Collections;

public class EigenOperations {
	
	
	
	public ArrayList Householder(Matrix A) {
		ArrayList QR = new ArrayList();
		
		Matrix A0 = A.copy();
		
		int col = A0.getColumnDimension();
		int row = A0.getRowDimension();
		
		ArrayList<Matrix> Householders = new ArrayList<Matrix>();
		
		for(int i = 0; i < row - 1; i ++) {
			for(int j = 0; j < col - 1; j++) {
				if(i == j) {
					// Make column vector and populate
					Matrix X = new Matrix(row - i, 1);
					double norm = 0;
					int currentRow = 0;
					for(int k = i; k < row; k++) {
						double a = A0.get(k, j);
						X.set(currentRow, 0, a);
						currentRow++;
					}
					norm = X.normF();
					Matrix Norm = new Matrix(row - i, 1);
					Norm.set(0, 0, norm);
					
					Matrix V = new Matrix(row - i, 1);
					V = X.plus(Norm);
					
					double vNorm = V.normF();
					Matrix VVt = V.times(V.transpose());
					Matrix vIdentity = A0.identity(row, row);
					
					for(int m = i; m < row; m++) {
						for(int n = j; n < row; n++) {
								vIdentity.set(m, n, VVt.get(m - i, n - j));
							}
						}
	
					Matrix Identity = A0.identity(row, row);
					Matrix H = Identity.minus((vIdentity).times(2/Math.pow(vNorm, 2)));
					Householders.add(H);
					A0 = H.times(A0);
				}
			}
		}
		
		Matrix R = A0;
		Matrix Q = Householders.get(0);
		for(int h = 1; h < Householders.size(); h++)
			Q = Q.times(Householders.get(h));
		QR.add(Q);
		QR.add(R);
		return QR;
		
	}
	
	public ArrayList getEigenvalues(Matrix A) {
		ArrayList<Matrix> QR = Householder(A);
		Matrix Q = QR.get(0);
		Matrix E = Q.times(A.times(Q.inverse()));
		
		ArrayList Eigenvalues = new ArrayList();
		
		for(int i = 0; i < E.getRowDimension(); i++) {
			Eigenvalues.add(E.get(i, i));
		}
		
		return Eigenvalues;
	}
	
	public ArrayList getEigenvectors(Matrix A) {
		ArrayList<Matrix> QR = Householder(A);
		Matrix Q = QR.get(0);
		
		ArrayList<Matrix> Eigenvectors = new ArrayList<Matrix>();
		for(int i = 0; i < Q.getRowDimension(); i++) {
			Matrix Eigenvector = Q.getMatrix(i, i, i, Q.getColumnDimension() - 1);

			Eigenvectors.add(Eigenvector);
		}
		
		return Eigenvectors;
	}
	
	public ArrayList SVD(Matrix A, int k) {
		ArrayList returnList = new ArrayList();
		
		Matrix AAt = A.times(A.transpose());
		Matrix AtA = A.transpose().times(A);
		
		ArrayList<Matrix> AAtQR = Householder(AAt); // V
		ArrayList<Matrix> AtAQR = Householder(AtA); // U
		
		Matrix AAtQ = AAtQR.get(0);
		Matrix AAtR = AAtQR.get(1);
		
		Matrix AtAQ = AtAQR.get(0);
		Matrix AtAR = AtAQR.get(1);
		
		ArrayList vList = getEigenvectors(AAtQ);
		ArrayList uList = getEigenvectors(AtAQ);
		
		Matrix V = new Matrix(((Matrix) vList.get(0)).getRowDimension(), vList.size()); // m x m
		for(int i = 0; i < vList.size() - 1; i++) {
			for(int j = 0; j < ((Matrix) vList.get(i)).getRowDimension() - 1; j++)
				V.set(i, j, ((Matrix) vList.get(i)).get(j, 0));
		}
		
		Matrix U = new Matrix(((Matrix) uList.get(0)).getRowDimension(), uList.size()); // n x n
		for(int m = 0; m < uList.size() - 1; m++) {
			for(int n = 0; n < ((Matrix) uList.get(m)).getRowDimension() - 1; n++) 
				U.set(m, n, ((Matrix) uList.get(m)).get(n, 0));
		}
			Matrix D = new Matrix(V.getColumnDimension(), U.getColumnDimension()); // m x n
			
			//Set D equal to the larger R
			Matrix matrixtoUse = AAt.copy();
			if(D.getRowDimension() > D.getColumnDimension())
				matrixtoUse = AtA.copy();
			else 
				matrixtoUse = AAt.copy();
			
			
				
			// Put SVDs in D from smaller R
			for(int x = 0; x < matrixtoUse.getRowDimension() - 1; x++) 
				for(int y = 0; y < matrixtoUse.getColumnDimension() - 1; y++) 
					if(x == y) 
						D.set(x, y, Math.sqrt(matrixtoUse.get(x, y)));
					
				
			ArrayList SVDs = new ArrayList();
			for(int g = 0; g < matrixtoUse.getRowDimension() - 1; g++)
				for(int h = 0; h < matrixtoUse.getColumnDimension() - 1; h++)
					if(g == h)
						SVDs.add(matrixtoUse.get(g,h));
			
			Collections.sort(SVDs);
			Collections.reverse(SVDs);
			
			ArrayList SVDcopy = new ArrayList();
			
			returnList.add(D);
			returnList.add(SVDs.get(0));
			if(k > SVDs.size()) 
				returnList.add((Integer) 0);
			else
				returnList.add(SVDs.get(k));
			
		return returnList;
	}
}
